^\ NIH Public Access 

i^^i Author Manuscript 

^'^HeP^ 7 Binmed Sci En^. Author manuscript; available in PMC 2014 March 24. 

Published in final edited form as: 
J Biomed Sci Eng. 2013 December ; 6(12A): 29-37. doi:10.4236/jbise.2013.612A005. 

Finite element modeling of the viscoelastic responses of the eye 
during microvolumetric changes 

Benjamin Cruz Perez^ Hugh J. Morris^ Richard T. Hart^ and Jun Liu^-^ * 

^Department of Biomedical Engineering, Tine OInio State University, Columbus, USA 

^Department of Ophthalmology, The Ohio State University, Columbus, USA 
Abstract 

A linear viscoelastic finite element model was built to investigate factors that influenced the 
intraocular pressure (lOP) elevations due to micro-volumetric changes in the eye at three different 
rates. The viscoelastic properties of the cornea and the sclera, including the instantaneous 
modulus, equilibrium modulus, and relaxation time constants, parametrically varied to examine 
their effects on lOP elevations at different rates of volumetric changes. The simulated responses 
were in good agreement with the previously reported experimental results obtained from porcine 
globes, showing the general trend of higher lOP elevations at faster rates. The simulations showed 
that all viscoelastic properties influenced the profile of the dynamic lOP due to volumetric 
changes, and the relative significance of a specific parameter was highly dependent on the rate of 
change. 
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1. INTRODUCTION 

Glaucoma is one of the leading causes of blindness worldwide [1]. It is known to be 
associated with elevated intraocular pressure (lOP), although the mechanisms of lOP- 
induced optic nerve damages are not well understood. The longstanding focus in the 
glaucoma field has been on the steady-state lOP [2], but it is also known that the lOP 
fluctuates during the day [3], as well as during normal physiological activities such as 
postural change [4,5], valsalva maneuver [6,7], water drinking [8], and eye movement [9]. 
Previous studies have shown that these short-term lOP fluctuations may be linked to 
glaucoma risk [8,10,11]. 

Because it has been very difficult to characterize the lOP-associated mechanical effects on 
ocular structures in vivo, computational modeling has been a helpful tool to describe and 
predict the biomechanical responses to lOP. Using finite element modeling of the posterior 
portion of the eye, Sigal et al. found that the material properties {i.e., the elastic modulus) of 
the peripapilary sclera played an important role in determining the deformation and stress 
states of the optic nerve head (ONH) [12], which was the primary site for early 
glaucomatous damages. 
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These findings have led to a growing interest in studying the biomechanical properties of the 
ocular shell and the association with glaucoma progression. Characterization of the 
mechanical properties of the cornea and sclera has been conducted, including uniaxial [13- 
16], biaxial [17], and inflation testing [18-22]. These studies showed that cornea and sclera 
exhibit the typical nonlinear, anisotropic, and viscoelastic mechanical behavior seen in many 
other biological tissues [23-27]. 

Very little has been done to understand the response of the eye to the short-term changes of 
volume and pressure. Recently, our group examined lOP elevations in porcine globes due to 
short-term micro-volumetric infusions. We found a strong rate-dependent response showing 
that a fast volumetric change induced a significantly higher lOP change in the same eye 
compared to slow volumetric changes [28]. The current study aims to build a finite element 
model of the whole corneoscleral shell (omitting the optic nerve head) to investigate how 
cornea and sclera viscoelastic properties influence the eye's response to microvolumetric 
changes and how each factor contributes to the rate-dependent response. 

The whole eye models have been developed in the past for understanding ocular physiology. 
Srodka [29] utilized an axisymmetric finite element of the whole eye model to evaluate the 
accuracy of the standard clinical method for lOP measurement, the Goldmann applanation 
tonometry. Anderson et al. [30] utilized a 3D hyperelastic whole eye model to investigate 
the influence of the boundary conditions at the cornea and sclera conjunction on the apical 
rise of the cornea during inflation. Our current study incorporated the general formulation of 
these previous whole eye models such as using an axisymmetric geometry and omitting the 
ONH. 

2. METHODS 

2.1. Finite Element lUlodeling 

An axisymmetric model of the corneoscleral shell of a porcine eye was constructed in 
COMSOL (v4.3, Burlington, MA, USA) using a series of ellipses assuming symmetry along 
the optic axis. The ellipsoidal building blocks allowed for realistic representation of the 
thickness variations in the corneoscleral shell. The central cornea thickness (CCT) was set to 
0.96 mm [31], with a radius of curvature of 7.5 mm, an assumed eccentricity of 0.5, and a 
white-to- white (WTW) diameter of 13.5 mm. The scleral radius of curvature was set to 12.8 
mm, the anterior thickness to 0.89 mm, the equatorial thickness to 0.58 mm, and the 
posterior thickness to 1.12 mm. Because of the oblate spheroidal shape of the porcine sclera, 
the ratio between the posterior radius and the equatorial radius was set to 0.86. These 
dimensions were the average measurements in porcine eyes (SiouxPreme Packing Co, Sioux 
City, lA) acquired for other studies in our laboratory. 

The cornea and sclera were divided by a radial line going through the center of the sclera 
ellipse. A small transition zone at the cornea-sclera junction was implemented with a finer 
mesh to avoid convergence difficultties at the boundary of material property change. A 
filling material was used to fill the inner space of the ocular shell and to couple the 
intraocular volumetric change to the mechanical responses of the shell (as described 
subsequently). Figure 1(a) shows the constructed geometry of the whole eye as described 
above, and Figure 1(b) illustrates the typical finite element mesh. 

The mesh was constructed using quadratic triangular elements. A displacement mesh 
convergence test was performed with five different mesh densities that monitored five 
different locations to ensure mesh accuracy. The mesh density chosen for this study had 
approximately 11,900 elements for a total of 54,000 degrees of freedom (Figure 1(b)). A 
time -dependent solver in COMSOL was utilized for solving the constitutive equations. 
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2.2. Simulated Infusions 



Our laboratory has previously reported an experimental infusion study on porcine whole 
globes at three different rates (the details of the experiments are provided in the referenced 
article) [28]. In the present study, the simulations were performed following the design of 
the infusion experiments. Specifically, three different infusion rates (as shown in Table 1) 
were implemented to represent the typical time scales seen in short-term lOP fluctuations 
(seconds to minutes) [28]. The total volume change was set to 15 )iL. 



2.3. Using Thermal Expansion to Simulate Intraocular Volumetric Change 



In order to simulate intraocular volume increments, it was computationally convenient to use 
a thermal expansion model although in reality there is no temperature dependence being 
studied. However, by assuming that the tissue and fluids inside the ocular shell are 
incomepressible, the volume change can be modeled as a simple uniform volumetric thermal 
expansion of an intraocular filling material, which can be implemented in the COMSOL 
software to conveniently provide mechanical coupling between the shell and the intraocular 
material. The thermal expansion model assumed that all mesh points within the enclosed 
filling material had a heat source and the total power (P) was computed by the following 
relation: 



where dV/dt is the volumetric change rate, Cp is the specific heat, pis the density and ay is 
the coefficient of thermal expansion for the filling material. The parameters for low density 
polyethylene (LDPE) was used to simulate the filling material because of its high coefficient 
of thermal expansion (ay= 4 x 10"^ 1/°C) [32]. The volumetric change rate (dV/df) resulted 
from the thermal expansion was set equal to the experimental infusion rate assuming no 
leakage. A fixed point was placed at the center of the filling material to ensure that it 
expands equally in all directions. The filling material was thermally insulated at the outer 
boundary to ensure conservation of energy. 



2.4. Constitutive Model 

The cornea and sclera were modeled using a nearly incompressible linear viscoelastic 
generaUzed Maxwell model. The stress-strain relationship is defined as follows [32]: 



where E^^i, is the equilibrium modulus, Ej is the relaxation modulus of the i-th branch, Zj is 
the time constant of the i-th branch, fis the strain, and sis the strain rate. A two-branch 
model {i.e., m = 2) was adopted [16]. The time constants (tj and t^) are denoted as the short 
term time constant ( z^) and the long term time constant ( t/) in the later discussions for clear 
differentiation. The instanttaneous {Einst) is defined as [14]: 



The branch relaxation moduli were set to be equal {i.e., E\ = E-i). Therefore, the branch 
relaxation moduli can be found as: 



(7 1 





(3) 



e 't=o 
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The instantaneous and equilibration moduli of the sclera were set to be 5 times of those of 
the cornea to ensure a stiffer response of the sclera commonly reported in experimental 
results and simulation studies [28,34,35]. All other parameters were allowed to vary 
independently. 

2.5. Parametric Study 

The values of the material parameters of ocular shell were varied within a range of the 
plausible values for the porcine eye (as shown in Table 2) to evaluate the effects that each 
parameter has on the short-term lOP elevations. To the best of our knowledge, experimental 
measurements of the viscoelastic properties of the porcine cornea and sclera including the 
instantaneous modulus, equilibrium modulus, and time constants have not been reported. 
The selected baseline values were estimated from approximating the model predictions with 
the average results obtained from the infusion experiments in porcine eyes [28]. 

3. RESULTS 

A comparison of the simulated lOP/infusion volume response (solid lines) and the average 
experimental responses in the porcine eyes (individual markers) is presented in Figure 2. All 
viscoelastic properties were set at the baseline values. Selected volume levels (11 levels 
evenly distributed between 0 to 15 |J,L) were used for plotting the experimental response for 
visual clarity. Most of the markers fell on the corresponding curves, suggesting a good 
agreement between the model simulation and the experimental results could be 
simultaneously achieved at all three rates. 

Figure 3 presents the effects of instantaneous modulus and equilibrium modulus on lOP 
elevations during micro-volumetric changes in the eye at different rates. The instantaneous 
moduli of the cornea and the sclera were varied proportionally (i.e., 1:5) while all other 
parameters were kept at their baseline. In general, lOP elevations were higher at higher 
modulus regardless of the infusion rates. Figure 3(a) shows that the lOP elevations increased 
almost linearly with the instantaneous moduli of the cornea and sclera at all rates. The slope 
however was different for different rates with a larger slope at the higher infusion rate. 
Figure 3(b) shows the effects of the equilibrium modulus. The response was less linear for 
the fast and intermediate infusion rates, particularly at the lower modulus range. 

Treating each data point in Figure 3(a) as an individual eye, the data points form a group of 
simulated eyes that only differ in instantaneous modulus. The relationships of the lOP 
elevations at different infusion rates in the group of eyes are plotted in Figure 4(a), showing 
that the lOP elevations at different infusion rates are predicted to be linearly correlated in a 
group of eyes that only differ in the instantaneous modulus. A similar plot for the eyes that 
only differ in equilibrium modulus (using data points in Figure 3(b)) is shown in Figure 
4(b). The relationships of lOP elevations at different infusion rates become nonlinear for this 
case. Figure 4(c) shows the experimental results obtained from 1 1 porcine eyes, showing 
that the lOP elevations at different infusion rates were approximately linearly correlated in 
the tested eyes (R = 0.99, p < 0.001, for fast and intermediate correlation; R = 0.89, p < 
0.001, for fast and slow correlation) [28]. 

Figure 5 presents the effect of the time constants on the lOP elevations due to micro- 
volumetric changes at different infusion rates. Figure 5(a) shows that the short term time 
constant affects lOP elevations during fast volumetric changes but not during intermediate 
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or slow volumetric changes. A larger short term time constant results in a larger lOP 
elevation for fast volumetric changes. Conversely, the long term time constant affects the 
intermediate and slow volumetric changes but not the fast (Figure 5(b)). A larger long term 
time constant results in larger lOP elevations for the slow and intermediate rates. 

Figure 6 shows the comparative influence of each of the corneal and viscoelastic properties 
on the lOP elevations, separately presented for each infusion rate. The x-axis shows the 
value for all parameters linearly scaled from its minimum value ( 1) to maximum value (+1). 
It is found that the instantaneous modulus has the largest influence for the fast and 
intermediate rates, but its influence diminishes for the slow infusion rate, for which the 
equilibrium modulus takes the dominant role. The time constants generally have smaller 
effects than the moduli, but their influence could be comparable to one of the two types of 
moduli depending on the rates. 

Figure 7 shows the distribution of the von Mises stresses and principal strains in the cornea 
and the sclera at the end of the fast infusion. The stresses and strains were apparently 
different in the cornea as compared to those in the sclera, due to the discrepancy in the 
material properties. The stresses and strains also showed variations through the thickness as 
contrasted to what a thin-shell model would predict [22]. 

4. DISCUSSION 

In this work, we have built a finite element model of the whole ocular shell to investigate 
factors influencing the short-term volume/pressure changes in the eye. We found that a 
simple linear viscoelastic model with two branches simulating the short-term and long-term 
relaxation processes could adequately reproduce the experimental data obtained from the 
infusion experiments performed in porcine eyes (Figm^e 2). Although the three different 
infusion rates used in this study are all within the time scale of short-term lOP fluctuations 
seen in the eye (i.e., a time scale shorter than a few minutes), there was a clear rate- 
dependent response suggesting that the eye responds differently to micro-volumetric 
changes that occur within seconds, lO's of seconds, or minutes. The larger the rate of change 
within this time scale, the higher the lOP elevations observed in both experiments and 
simulations. This rate-dependence confirmed the viscoelastic responses of the ocular shell, 
as well as the relevant range for the time constants of the relaxation processes in the ocular 
shell (i.e., seconds to minutes), which were found using stress relaxation tests on tissue 
strips [16]. Both experiments and simulations reported a linear AlOP/volume relationship, 
which was also experimentally observed by Pierscionek et al. in fresh porcine eye [36]. 

Based on the finite element model, a higher instantaneous modulus, i.e., a stiffer response to 
a step load, would result in a larger final AIOP for a given infusion volume regardless of the 
rate of change (Figure 3(a)). A similar positive relationship between AIOP and equilibrium 
modulus was also shown (Figure 3(b)). The relationship between instantaneous modulus and 
AIOP relationship was linear (i.e., a constant slope) while there was a degree of nonlinearity 
in the relationship between equilibrium modulus and AIOP, although that was mostly seen in 
the fast rate response. 

In our previous experimental study, we found a strong correlation between the lOP 
elevations measured at different rates across different eyes (Figure 4(c)), that is, an eye that 
had a larger AIOP at the fast rate would also have a larger AIOP at the slower rates. 
Interestingly, this correlation was strongly linear (R = 0.99 or 0.89) according to the 
experimental data. This linear correlation could be better predicted by the variance in 
instantaneous modulus (Figure 4(a)) rather than the variance in the equilibrium modulus 
which corresponded to a more nonlinear response (Figure 4(b)). 
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The short term time constant had a positive relationship with AIOP during fast infusions, but 
a flat response for intermediate and slow infusions (Figure 5(a)). This time constant 
represents the time scale of the fast relaxation processes in the tissue. Because of the 
selected range of this constant in the present study (0.05 to 0.65 sec), it is expected that it 
does not influence the intermediate and slow rates whose time scales were far beyond this 
range. The long term time constant represents the slower relaxation process. As shown in 
Figure 5(b), the AIOP in the intermediate and slow infusions were noticeably affected by 
this time constant, but the AIOP in the fast infusion had a flat response because its short 
duration. 

As discussed above, all viscoelastic properties affected the volume/pressure response, and 
the relative significance of each specific parameter was rate-dependent. For example, at the 
faster rates (e.g., the fast and intermediate infusion rates used in the present study), the 
instanttaneous modulus had the largest influence on the magnitude of AIOP (Figures 6(a) 
and (b)). As the rate decreased, the equilibrium modulus started to dominate (Figure 6(c)). In 
addition, the short-term time constant had minimal effect on intermediate and slow 
infusions, but significant effect on fast infusion; while the long-term time constant showed 
the opposite trend. 

Previous computational work has shown that the posterior sclera is a major biomechanical 
structure that strongly influences the mechanical states of the optic nerve head where the 
glaucomatous injuries occur [12]. The results from the present study showed that the 
viscoelastic properties of the corneoscleral shell could also influence the dynamic profile of 
lOP. For example, a larger instantaneous modulus or equilibrium modulus (i.e., a stiffer 
corneoscleral shell) would result in higher lOP fluctuations due to micro-volumetric changes 
in the eye during daily activities. Although the pathophysiolgical consequences of lOP 
fluctuations are not fully understood, it is generally believed that they are detrimental to the 
delicate tissue structures in the eye, especially when they are at high magnitudes and fast 
acting [37,38]. It is therefore important to understand and characterize the viscoelastic 
properties of the ocular shell and their relationship to the dynamic lOP. 

This work has several limitations. First, the eye is actually a "leaky" shell with a pressure- 
dependent outflow [39]. For simplicity, this pressure-dependent outflow was not 
implemented in the current model. Instead, a net flow rate was used as the infusion rate. This 
could have introduced certain inaccuracy in the calculations of the AIOP; however, this error 
was likely small due to the short time duration of the simulated experiments. Second, the 
relaxation moduli for the branches in the Maxwell model were set to the same values. 
Although this assumption has been used in the past [14], future experimental studies as well 
as modeling work are needed to verify this assumption and investigate its implications. 
Third, a linear viscoelastic model was used and the reported nonlinearity in the stress-strain 
relationship of the cornea and sclera was neglected. Because the corneoscleral shell is 
primarily a collagenous tissue, mechanical nonlinearity is expected and has been confirmed 
in the past. Constitutive modeling that accounts for the nonlinear responses would thus be 
more desirable [40^4]. However, in the present study, the lOP was always above the 
baseline of 15 mmHg during the simulated infusion studies, which likely positioned the eye 
close to or beyond the toe region of the corneoscleral response [28]. Therefore, the linear 
approximation, which simplifies the model, is likely sufficient in the present study. Fourth, 
previous studies have shown the importance of collagen microstructure (for example, 
anisotropic collagen fiber alignment) in affecting the stress and strain distributions in the 
cornea and sclera [20,45^8]. Future work should investtigate how heterogeneous 
microstructural and material properties could affect the eye's response to dynamic micro- 
volumetric changes. 
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5. CONCLUSION 

In conclusion, this study utilized finite element models of a porcine corneoscleral shell to 
investigate the effects of the shell's viscoelastic properties on the lOP elevations due to 
volumetric changes in the eye. The major finding was that the instantaneous and equilibrium 
moduli, as well as the short-term and long-term time constants, had an effect on lOP 
elevations, and the relative significance of the specific parameters was highly dependent on 
the rate of change. 
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Figure 2. 

The lOP elevations (AIOP) as a function of the infused volume as predicted by the model 
(solid hnes) or obtained from the experiments (markers). 
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Figure 3. 

The effects of corneal and scleral biomechanical properties on lOP elevations at different 
infusion rates, (a) Effects of instantaneous modulus, (b) Effects of equilibrium modulus. 
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Figure 4. 

The relationships of lOP elevations at different infusion rates in: (a) A group of simulated 
eyes that only differ in instantaneous modulus; (b) A group of simulated eyes that only differ 
in equilibrium modulus; and (c) 1 1 porcine eyes used in the experiments [38]. 



J Biomed Sci Eng. Author manuscript; available in PMC 20 14 March 24. 



Perez et al. 



Page 14 




7 H 1 1 1 1 

0.05 0.2 0.35 0.5 0.65 
Short Term Time Constant (sec) 

(a) 



— Fast — Intermediate — Slow 
14.5 H 




5.5 T 1 1 1 1 

18 43 68 93 118 

Long Term Time Constant (sec) 

(b) 



Figure 5. 

The effects of (a) the short term time constant and (b) the long term time constant on lOP 
elevations at different infusion rates. 
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Figure 6. 

Effects of all corneal scleral factors on lOP elevations due to micro-volumetric changes at 
the (a) fast, (b) intermediate, and (c) slow rate. 
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Figure 7. 

Color maps of von Mises stress (a) and first principal strain (b) in the cornea and the sclera 
at the end of the fast infusion (15 |aL/s for 1 second). 
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Table 1 

The infusion rates and total duration used in the experiments and simulations. 





Fast 


Intermediate 


Slow 


Rate 


15 (iL/s 


1 [iL/s 


0.1 |iL/s 


Duration 


1 s 


15 s 


150 s 
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Table 2 

Material parameters used in the model simulations. 

Sclera Cornea 

Material properties 

Baseline Range Baseline Range 

Instantaneous modulus 5.5 MPa 3.1 - 7.9 MPa 1.1 MPa 0.62 - 1.6 MPa 

Equilibrium modulus 1.325 MPa 0.325 - 2.325 MPa 0.265 MPa 0.065 - 0.465 MPa 

Short term time constant 0.35 sec 0.05 - 0.65 sec 0.35 sec 0.05 - 0.65 sec 

Long term time constant 68 sec 18 -128 sec 68 sec 18 -128 sec 
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